function [ out ] = fun_sigma_est(xx,theta,combos_all_unique,m_regr,J)

temp = zeros(1,size(combos_all_unique,1));

for i=1:size(combos_all_unique,1)
     
    ind_other = perms(1:J-1);
        
    combos_x_other = combos_all_unique(i,2:J);
    combos_x_other_r = combos_x_other(ind_other);
    
    combos_z_other = combos_all_unique(i,J+2:2*J);
    combos_z_other_r = combos_z_other(ind_other);
    
    v = ([repmat(combos_all_unique(i,1),[size(ind_other,1),1])   combos_x_other_r...
                repmat(combos_all_unique(i,J+1),[size(ind_other,1),1]) combos_z_other_r]) - 1;
   
    v = unique(v,'rows');       % KEY to avoid counting multiple times the combos that correspond to only one of combos_all        
            
    temp(i) = double(sum(prod(BernsteinElement_2(v,m_regr,xx'),3)',2));

end

out = temp*theta;

end